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SUMMARY 


An  on-line  technique  is  presented  for  the  identification  of  rotor  blade 
modal  damping  and  frequency  from  rotorcraft  random  response  test  data.  The 
identif ication  technique  is  based  upon  a  recursive  maximum  likelihood  (RML) 
algorithm,  which  is  demonstrated  to  have  excellent  convergence  characteristics 
in  the  presence  of  random  measurement  noise  and  random  excitation.  The  RML 
technique  requires  virtually  no  user  interaction,  provides  accurate  confidence 
bands  on  the  parameter  estimates,  and  can  be  used  for  continuous  monitoring  of 
modal  damping  and  frequency  during  wind  tunnel  or  flight  testing. 

Results  are  presented  from  simulation  random  response  data  which  quantify 
the  identified  parameter  convergence  behavior  for  various  levels  of  random  excita¬ 
tion.  The  data  length  required  for  acceptable  parameter  accuracy  is  shown  to 
depend  upon  the  amplitude  of  random  response  and  the  modal  damping  level.  Random 
response  amplitudes  of  1.25  degrees  to  .05  degrees  are  investigated.  The  RML 
technique  was  applied  to  hingeless  rotor  test  data  from  the  NASA  Langley  Research 
Center  Helicopter  Hover  Facility.  The  inplane  lag  regressing  mode  damping  was 
correctly  identified  at  different  rotor  speeds.  The  identification  from  the  test 
data  agreed  with  the  simulation  results  and  with  other  available  estimates  of 
frequency  and  damping. 


INTRODUCTION 

The  determination  of  rotor  blade  modal  damping  for  lightly  damped  rotors 
during  wind  tunnel  or  flight  testing  is  an  important  safety  requirement.  Analytic 
methods  to  predict  stability  are  often  employed  before  testing,  however  the 
accuracy  of  analytic  methods  are  not  always  reliable,  particularly  for  new  rotor 
designs  or  configuration  modifications.  Estimation  of  rotor  blade  mode  damping 
from  test  data  is  thus  required  during  wind  tunnel  testing  where  the  possibility 
of  low  damping  levels  may  occur.  Current  methods  for  estimation  of  modal  damping 
from  test  data  are  based  upon  use  of  the  Fast  Fourier  Transform  (FFT) .  This 
technique  is  often  referred  to  as  the  moving  block  FFT  method.  This  method  has 
been  applied  to  rotor  data  obtained  either  from  an  initial  condition  decay 
response  obtained  after  removal  of  a  forced  excitation,  or  from  random  response 
data  (ref.  1,  ref.  2  and  ref.  3).  The  random  response  data  is  first  transformed 
to  an  initial  condition  decay  response  by  calculating  an  ensemble  average  of 
segments  of  the  random  response  time  history.  Application  of  the  moving  block 
FFT  approach  is  then  performed.  This  technique  has  been  referred  to  as  the  moving- 
block/randomdec  method  (ref.  1). 

During  wind  tunnel  testing,  it  would  be  desirable  if  damping  and  frequency 
could  be  computed  without  the  use  of  a  forced  input  excitation  and  be  computed 
in  real  time  along  with  confidence  bands.  Although  the  FFT  technique  has 
received  widespread  acceptance  by  the  rotorcraft  community,  practical  usage  of 
this  technique  requires  considerable  user  interaction  and  is  not  a  real  time  method. 
Estimation  of  damping  and  frequency  requires  several  minutes  and  continued  user 


interaction.  In  addition,  wind  tunnel  applications  are  usually  done  with  a  forced 
input  excitation  to  obtain  the  decay  response  time  history. 


The  RMI.  technique  presented  in  this  report  requires  virtually  no  user  interac¬ 
tion,  and  provides  on-line  continuous  estimates  of  damping  and  frequency  from 
random  response  data.  In  addition,  accurate  confidence  bands  (+  la  bands)  are 
also  recursively  computed  which  indicate  the  required  data  length  for  accurate 
parameter  estimation.  The  technique  can  be  used  on  very  small  amplitude  random 
response  data  which  often  is  the  case  in  wind  tunnel  tests.  During  wind  tunnel 
testing  it  is  common  to  have  data  continuously  generated  while  steady  trim  condi¬ 
tions  are  achieved,  and  this  data  is  suitable  for  parameter  identification. 

Thus,  the  identif ication  procedure  need  not  interfere  with  other  wind  tunnel  test 
objectives.  Although  the  RML  technique  presented  is  developed  for  any  number  of  modes, 
including  closely  spaced  or  overlapping  modes,  the  simulation  and  test  results 
presented  in  this  report  focus  on  the  single  mode  case. 

Helicopters  can  encounter  many  different  types  of  rotor  dynamic  instability 
during  wind  tunnel  testing.  Instability  can  arise  from  pure  mechanical  coupling 
as  in  ground  resonance  or  from  aerodynamic  coupling  in  forward  flight.  Rotor 
blades  are  designed  such  that  the  rotor  blade  mode  frequencies  lie  between 
multiples  of  the  rotor  normal  operating  RPM.  This  fact  places  the  modes  of 
interest  between  a  known  band  of  frequencies  and  allows  isolation  of  these  modes 
with  the  use  of  a  bandpass  filter.  (This  is  in  contrast  to  fixed  wing  flutter 
where  the  modal  frequencies  are  not  as  well  defined.)  Prior  knowledge  of  the 
approximate  frequency  band  allows  for  pretest  seleccion  of  the  bandpass 
filter  passband,  which  in  turn  permits  automated  use  of  the  identification 
algorithm.  The  identification  method  identifies  both  damping  and  frequency  of 
the  modes  within  the  preselected  passband. 

Recursive  maximum  likelihood  methods  for  the  identification  of  parameters 
from  input/output  data  have  been  extensively  researched  (ref.  4  and  ref.  5). 

However,  the  application  to  real-time  identification  of  rotor  blade  damping  has 
not  been  previously  done.  These  real-time  methods  are  based  upon  an  input /output 
discrete  transfer  relationship  or  autoregressive  moving  average  (ARMA)  representa¬ 
tion.  This  representation  leads  to  a  simple,  accurate,  and  numerically  efficient 
identification  algorithm.  Other  methods  for  real-time  identification  of  para¬ 
meters  such  as  least  squares  or  instrumental  variable  produce  either  biased 
estimates  in  the  presence  of  measurement  noise  or  require  extensive  user  pre¬ 
tuning  for  proper  convergence.  The  recursive  maximum  likelihood  procedure  selected 
for  the  research  reported  herein  has  excellent  convergence  properties  and  requires 
no  user  pretuning.  This  technique  is  referred  to  as  the  RML2  method  in  reference 
4  and  reference  5. 

Flutter  boundary  prediction  has  been  performed  in  reference  6  for  a  fixed 
wing  aircraft  using  nonrecursive  maximum  likelihood  identification  in  conjuction 
with  Jury's  determinate  method  to  predict  the  stability  boundary.  The  present 
research  differs  from  the  techniques  used  in  reference  6  in  that  recursive  on¬ 
line  maximum  likelihood  procedures  are  used.  In  addition,  the  objective  in  this 
research  is  to  provide  continuous  on-line  estimation  of  modal  damping  as  wind 
tunnel  testing  is  performed. 

The  real  time  identification  procedure  is  described  in  this  report.  An 
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investigation  is  conducted  using  simulation  data  to  determine  the  convergence 
characteristics  of  the  algorithm  and  identification  accuracy  with  random  input 
excitation  and  measurement  noise.  The  effect  of  bandpass  filter  bandwidth  on 
identification  accuracy  is  quantified.  A  lightly  damped  blade  mode  with  2% 
critical  damping  and  moderately  damped  mode  of  10%  critical  damping  is  used  in 
the  investigation.  The  algorithm  is  applied  to  hingeless  rotor  random  response 
test  data  obtained  from  the  NASA  Langley  Research  Center  Hover  Test  Facility. 
Test  data  is  analyzed  using  the  algorithm  and  damping  and  frequency  estimates 
are  compared  with  other  available  predictions. 
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SYMBOLS 


+  bnz‘ 


coefficients  of  the  characteristic  equation  for  a  dynamic  system  of 
order  n,  j  =  1,  n,  also  defined  in  equation  2.9 
stability  matrix  appearing  in  linear  state  equation  2.5 

polynomial  in  the  shift  operator  z~^ ,  A(z)  =  1  +  a^z--*-  +  a2Z~^  + - +  anz 

control  matrix  appearing  in  linear  state  equation  2.5 

coefficients  associated  with  the  control  input  in  the  ARMA  representation 

polynomial  in  the  shift  operator  z~^ ,  B(z)  =  b^z--*-  +  b2Z-^  + - 

cyclic  control  input  to  rotor,  deg. 
filter  bandwidth  used  in  the  digitial  filter,  Hz 

coefficients  associated  with  the  innovations  in  the  ARMA  representation 

polynomial  on  the  shift  operator  z~^ ,  C (z)  =  1  +  c^z-^  +  C2z'^  +  -  +  c 

covariance  of  the  parameter  estimate 

covariance  of  the  parameter  estimate  corrected  for  nonwhite  innovations 
natural  logrithm  basis,  e  =  2.71828 

generalized  force  scaled  by  generalized  inertia  for  the  ith  mode 
control  matrix  appearing  in  the  discrete  state  equation  2.7 
measurement  matrix  appearing  in  the  linear  state  equation  2.7 
weights  used  in  the  digital  filter 
discrete  data  point  index,  also  mode  number  index 


complex  variable  j  --  /-I,  also  data  sample  shift  index 
Kalman  filter  gain 

gain  matrix  used  in  the  recursive  maximum  likelihood  algorithm 
negative  log  likelihood  function 
likelihood  function 

number  of  data  points  on  either  side  of  current  point  fitered  by  the 
digital  filter 

vector  of  mode  shapes  at  a  radial  station 
number  of  data  points 


-n 


number  of  coefficients  in  ARMA  model  for  a^,  b^,  c^;  i  =  l,n 
parameter  covariance  matrix  used  in  the  recursive  maximum  likelihood 
algorithm 

initial  parameter  covariance 
nondimens ional  rotor  frequency 

generalized  rotor  coordinate  associated  with  ith  mode 
root  mean  square 
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Laplace  transform  variable  s  -  a  +  jco 
time  variable,  sec,  also  discrete  time  index 
data  sample  time  interval,  sec 
control  input 

discrete  white  gaussian  random  measurement  noise 
loss  function  minimized  in  the  RML  algorithm 
discrete  white  gaussian  random  process  noise 
zero-mean  continuous  white  gaussian  process  noise 
state  vector 

vector  of  sensitivity  equations  used  in  the  RML  algorithm 
measurement  variable  (rr  vector) 
digital  filtered  data 

z-transform  shift  variable,  also  modal  response  variable 

real  and  imaginary  parts  of  the  complex  variable  Z,  2  =  +  j 

innovations 

process  noise  gain  matrix 
innovation  covariance 
convergence  parameter 
innovations 

discrete  state  transition  matrix 
parameter  vector 
initial  parameter  vector 

standard  deviation,  also  used  as  real  part  of  eigenvalue 

standard  deviation  of  process  noise  and  measurement  noise  respectively 

standard  deviation  of  damping  parameter  estimate 

frequency,  rad/ sec. 

undamped  natural  frequency,  rad/sec 

damped  frequency  associated  with  inplane  mode,  rad /sec 

cucoff  and  termi'.iation  frequency  used  in  the  digital  filter,  rad/sec 

rotor  speed,  revolutions  per  min 

modal  damping,  %  critical 

denotes  transpose 

denotes  time  derivative 

denotes  estimated  quantity 

partial  derivative  with  respect  to  0 
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THEORETICAL  BACKGROUND 

This  section  presents  the  problem  description  and  overall  approach  used  in 
the  identification  of  blade  modal  damping  and  frequency  from  test  data.  The 
real-time  maximum  likelihood  identification  algorithm  is  then  presented  followed 
by  the  bandpass  filter  used  to  separate  the  desired  modes  from  the  measured  test 
data. 


Problem  Description  and  Overall  Approach 

The  real  time  identification  of  helicopter  blade  damping  and  frequency  is 
based  upon  the  blade  modal  response  to  random  excitation.  This  response  is 
assumed  to  be  modeled  by  a  set  of  coupled  second  order  equations  with  a  single 
output  measurement  occurring  at  discrete  time  increments.  The  rotor  blade  response 
contains  frequencies  occurring  at  multiples  of  the  rotor  frequency  (lp,  2p,  etc.) 
and  the  various  blade  modes.  Bandpass  filtering  is  employed  to  separate  the  main 
modes  of  interest  from  the  rotor  lp,  2p,  etc.  frequencies  and  other  modes.  The 
resulting  response  data  can  then  be  considered  to  be  a  low  order  model  of  N 
coupled  second  order  equations.  Rotor  blade  response  in  this  passband  will 
typically  consist  of  froti  one  to  three  blade  modes. 


The  general  response  of  the  ith  blade  mode  of  the  rotor  in  the  rotating 
system  can  be  represented  as 


q.  +  2£. w  q.  +w  q  =  F. 
ni  1  nt  Mi  n^  i  1 


(2.1) 


where  q^  is  the  generalized  coordinate  and  is  the  generalized  force  (scaled 
by  the  generalized  inertia)  of  the  ith  mode.  The  vertical  deflection,  Z(t)  at  a 
radial  station  is  the  sum  of  the  modes 


#  modes 

z(t)  =  E  ^ 

i=l  1 


(2.2) 


where  is  a  vector  of  mode  shapes  at  station  r.  The  variable  Z(t)  can  also  be 
used  to  represent  the  blade  bending  moment  at  station  r  as  measured  by  strain 
gauage  instrumentation.  Defining  q^  =  M^q^  and  substituting  into  Eq.  (2.1)  and 
Eq.  (2.2)  eliminates  the  mode  shape  vector  from  the  vertical  deflection  Z. 
Equations  (1.1)  and  (2.2)  thus  becomes 

’i  +  2h\  ’i  +  “n  ‘  'Fi  <2'3> 

and 

it  modes 

z(t)  =  E  h  <2-4) 

i=l  1 

where  is  scaled  by  the  corresponding  mode  shape  M^. 

A  canonical  representation  for  the  deterministic  single  input-single  output 


5 


] 


w 


system  described  by  Eq.  (2.3)  and  Eq.  (2.4)  (ref.  7)  can  be  described  in  state 
form  as 


x  =  A  x  +  B’u 


(2.5) 


where 


(2.6) 


B' 


and  x  is  a  n  x  1  state  vector  of  normalized  coordinates  and  velocities,  u  is 
the  single  input  and  y  is  the  single  output. 

In  the  presence  of  random  process  and  measurement  noise  the  continuous  system 
Eqs .  (2.5)  and  (2.6)  can  be  represented  by  the  discrete  stochastic  system  of  linear 
difference  equations 


x(i+l)  =  <j>x(i)  +  Gu(i)  +  Tw(i) 


(2.7) 


y(i)  =  Hx(i)  +  v(i) 

where  i  is  the  index  for  the  discrete  sample  time.  The  discrete  time  process 
noise  w(i)  and  measurement  noise  v(i)  are  assumed  to  be  zero  mean  white  sequences 
with  covariance  matrices  Q  and  R  respectively. 

By  use  of  the  Kalman  Filtering  theorem  (ref.  8  and  ref.  9)  the  stochastic 
model  of  Eq.  (2.7)  can  be  represented  by  an  equivalent  form  as 


i 


i 


! 
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x(i-f-l)  =  <}>  x(i)  +  Gu(i)  4-  Kv(i) 
y(i)  =  H  x(i)  +  v(i) 


(2.8) 


where  x(i)  denotes  the  conditional  mean  of  x(i)  Riven  measurements  y(i-l),  y(i-2), 
...»  and  (v(i)}  is  a  sequence  of  independent  equally  distributed  random  variables 
with  zero  mean  and  covariance  X2.  The  representation  of  Eq.  (2.8)  is  an  innova¬ 
tion  representation  (ref.  8  and  ref.  9)  of  the  process  and  the  sequence  (v(t)} 
is  called  the  innovations.  In  steady  state  the  filter  gain  matrix  K  is  constant 
and  combined  with  the  innovations  covariance  X^,  replaces  the  noise  statistics  Q 
and  R  of  Eq.  (2.7).  The  innovations  sequence  v(i)  replaces  both  the  process  noise 
and  measurement  noise. 


The  innovation  model  of  Eq.  (2.8)  takes  on  a  particularly  simple  form  for  a 
single  input-single  output  (SISO)  model.  The  discrete  SISO  form  or  autoregressive 
moving  average  (ARMA)  representation  is  obtained  using  the  shift  operator  q_i,  where 
q~1x(i+l)  =  x(i) .  Application  of  the  shift  operator  to  Eq.  (2.8)  yields  the  SISO 
ARMA  representation 


n-1  n 

2^  aiy(1_j)  =  V  {b  u(j-l)  +  c  v(i-j)}  (2.9) 

j=0  J  H  2  2 


where  a0=l  and  c0  =  1,  n  is  the  order  of  the  autoregressive  part,  and  also  the 
order  of  the  moving  average  part.  In  this  description  there  is  no  concept  of 
state  and  a  n  state  vector  produces  n  lags  y(i-j),  j=0,  1,  2,... n-1. 


The  ARMA  representation  of  Eq.  (2.9)  provides  a  complete  description  of  the 
SISO  stochastic  process  of  Eq .  (2.7).  The  unknown  parameters  of  the  model  in 
Eq.  (2.9)  are  the  parameter  aj ,  bj ,  Cj  and  the  scalar  innovation  covariance.  The 
unknown  parameter  vector  0  is 


0  =  (a^,  a2, 


V  c2’ 


c 

n 


(2.10) 


n'  i'  z'  n' 

along  with  the  scalar  innovation  covariance  denoted  by  X^.  The  left-hand  side  of 
Eq.  (2.9)  contains  the  coefficients  of  the  characteristic  equation. 


The  damping  and  frequency  of  the  n-modes  are  obtained  from  the  roots  of  the 
characteristic  equation. 


„n-l 

a^  Z  + 


„n-2 
Z  + 


=  0 


(2.11) 


where  Z  is  the  discrete  Z  transform  variable. 


The  general  ARMA  model  given  by  Eq .  (2.9)  can  thus  be  used  to  identify  the 
coefficients  aj,  j*>l,  n-1  of  a  general  nta  order  system  from  the  single  measure¬ 
ment  y.  The  stability  of  the  system  is  determined  from  the  roots  of  the  discrete 
characteristic  equation.  To  obtain  the  continuous  damping  and  frequency  the 
mapping  from  the  discrete  z-transform  domain  to  the  continuous  s-transform  Laplace 
transform  domain  is  then  required. 
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The  ARMA  representation  of  Eq.  (2.9)  is  more  suitable  for  on-line  identifica¬ 
tion  than  the  complete  state  space  model  of  Eq .  (2.7)  or  the  state  innovations 
model  of  Eq.  (2.8).  The  ARMA  representation  has  fewer  parameters  and  the  noise 
statistics  are  completely  characterized  in  the  moving  average  part  of  Eq.  (2.9) 

(i.e.  the  c  coefficients)  and  innovation  covariance  A^.  Thus  in  order  to  identify 
the  coefficients  of  the  characteristic  equation  (i.e.  the  a  coefficients),  it  is 
also  required  to  identify  the  c  coefficients  and  innovation  covariance  A^.  This 
then  provides  a  complete  description  of  the  stochastic  process,  from  which  the 
modal  damping  and  frequency  can  be  determined.  Also,  if  a  known  control  input  is 
applied  the  b  coefficients  of  Eq.  (2.9)  must  also  be  identified. 

The  overall  identification  procedure  is  shown  in  block  diagram  form  in  figure 
1.  The  data  y  (and  control  input  u(i)  if  available)  is  discretized  and  bandpass 
filtered  to  obtain  the  desired  rotor  blade  modes.  The  filtered  data  is  passed 
through  the  recursive  maximum  likelihood  identification  algorithm  (discussed  in 
the  next  section)  and  the  identified  parameters  a^  represent  the  coefficients  of 
the  discrete  characteristic  equation.  Discrete  roots  (eigenvalues)  are  computed 
from  the  characteristic  equation.  The  modal  damping  and  frequency  are  then  obtained 
by  transforming  the  discrete  roots  to  the  continuous  Laplace  transform  domain,  ""he 
RML  algorithm  computes  both  the  parameters  and  standard  deviations  of  the  para  cr 
estimates  as  noted  in  figure  1. 

The  discrete  roots  are  obtained  from  the  characteristic  equation  by  use  o  iy 
standard  root  finding  algorithm  of  a  polynomial  equation.  Transformation  to  t 
continuous  domain  is  accomplished  by  defining  the  discrete  complex  root  Z  and 
continuous  root  s 


z  -  ZR  +  j  Zx 

(2.12) 

s  =  a  +  j  a) 

(2.13) 

and  using  the  discrete  to  continuous  transformation 

eST  =  Z 

(2.14) 

where  T  represents  the  discrete  sample  time  interval. 

Substituting  Eq.  (2.12)  and  (2.13)  into  Eq .  (2.14)  yields 


Expressing  Eq.  (2.15) 
hand  sides  yields 


uT  +  jwT  „  ,  „ 

« 


(2.15) 


in  terms  of  amplitude  and  phase  and  equating  right  and  left 


e2<TT  =  Z2  +  Z2 

-1  ZI 
u)T  =tan  — - 

'R 


(2.16) 

(2.17) 


Taking  the  natural  logrithm  of  Eq.  (2.16)  and  solving  for  a  and  m  in  Eq .  (2.16)  and 
(2.17)  yields  the  continuous  damping  and  frequency 


8 


(2.18) 


u>  =  y  cos_1(-a1/2v/a2)  (2.25) 

The  next  section  describes  the  details  of  the  recursive  maximum  likelihood 
identification  algorithm  for  computing  the  coefficients  of  the  ARMA  represention 
of  Eq.  (2.9). 

Recursive  Maximum  Likelihood  Identification  Algorithm 

The  origin  of  this  method  is  based  upon  the  off-line  maximum  likelihood  method 
of  Astrom-Bohlin  (ref.  10).  The  recursive  version  is  described  in  Soderscrom  (ref. 
11).  This  method  is  selected  for  the  rotorcraft  blade  mode  identification  due  to 
its  excellent  convergence  properties.  Reference  4  and  Reference  5  describes  the 
recursive  version  in  comparison  to  other  methods. 

The  off-line  maximum  likelihood  version  attempts  to  minimize  the  negative  log 
likelihood  function 

1  N  2 

L'  =  -log  L(0,X)  =  — -j  V  e  +  N/2  log  X  +  N/2  log  2tt  (2.26) 

2X  t=l 

where 
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and  A ^  represents  the  innovation  covariance.  The  ARMA  model  of  Eq.  (2.9)  is 
rewritten  using  the  operator  C(z) ,  A(z)  and  B(z) 

C(z)  et  =  A(z)  yt  -  B(z)  ut  ,  t  =  1,2, ...N  (2.28) 


where  the  innovations  is  denoted  by  and  t  denotes  the  discrete  sample  time. 

Since  the  prediction  error  Eq.  (2.28)  is  nonlinear  in  the  parameters  c,  a 
iterative  solution  is  required  to  minimize  the  likelihood  function  L'. 

2 

Optimization  of  L'  with  respect  to  A  yields 


1/N 


(2.29) 


Therefore  maximum  likelihood  identification  of  0  is  equivalent  to  minimization 
of 


v(0)  =  £  e2t 

t=l 

subject  to  the  prediction  error  model 

£t  =  [C(z)]'1  [A(z)  yt  -  B(z)  ut] 

An  approximate  solution  is  obtained  by  linearization  of  c 


(2.30) 


(2.31) 


et  =  e<0)  +  ||  (6)  (6  -  0) 


(2.32) 


where  the  linearization  is  performed  about  the  estimate  9. 

Eq.  (2.32)  is  a  linear  prediction  error  equation  and  therefore  the  minimization 
of  Eq.  (2.30)  subject  to  Eq.  (2.32)  can  be  obtained  using  the  recursive  least  square 
solution  (ref.  12) 

Therefore  the  recursive  solution  is 


A 


where  e 


t+1 


st+l  ■  K  *  Kt+1  ?t+l 

(2.33) 

•w  ■  pt  Vi/(xc  +  W  pt  W 

(2.34) 

pt+1  ■  [pt  -  Vi  ?+!  pt1/x« 

(2.35) 

is  computed  from 
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-c1  1  .  .  .  0 

1 

7” 

“b7 

O 

rH 

O 

CM 

<U 

1 

0 

x(t+l)  = 

•  •  • 

x(t)  + 

yt+i + 

yt  ~ 

:  1 

-c  0 

0 

/s 

a 

A 

b 

n 

..  „ 

L 

ELI 

L 

?J 

(2.36) 


where  x1 (t)  =  £(t) . 

1  3e  * 

The  sensitivity  xt+2  =  '30 -  is  computed  according  to 


(2.37) 


The  factor  Ac  in  Eq.  (2.34)  and  (2.35)  is  used  to  improve  convergence  rate 
(ref.  4)  and  is  obtained  from 

A*;+1  -  C  +  1  -  C  (2.38) 

COCO 

where  X°  =  .95  and  C  =  .99.  The  values  of  X°  and  C  are  the  same  as  those  used  in 

CO  CO 

reference  4  which  were  found  to  yield  good  convergence. 


11 


1 


To  track  slowly  time  varying  parameters  the  convergence  factor  Ac  is  set  to 
a  constant  value  less  than  unity  for  all  time.  This  is  equivalent  to  exponential 
weighting  of  the  least  square  criterion  Eq.  (2.30)  and  is  fully  discussed  in 
reference  4. 

Eqs.  (2.33)  through  Eq.  (2.38)  constitutes  the  recursive  maximum  likelihood 
algorithm. 


Bandpass  Filter 

The  selection  of  the  desired  modes  for  identification  is  accomplished  with  a 
zero  phase  shift  digital  filter,  due  to  Graham  (ref.  13).  This  filter  is  not 
recursive  but  uses  numerical  convolution  of  the  raw  data  to  obtain  the  filtered 
data.  This  filter  can  be  used  in  near-real  time  by  shifting  the  filter  weight 
function  hj  one  data  point  ahead  as  each  new  raw  data  point  is  received.  The 
actual  filtered  output  requires  a  finite  segment  of  data  ahead  in  time  of  the 
current  data  point  being  filtered.  This  results  in  a  slight  time  delay 
before  the  estimated  coefficients  are  available.  Using  a  sample  rate  of  .01 
seconds  this  represents  a  delay  of  between  1  and  2  seconds  before  the  filtered 
output  is  obtained.  A  recursive  bandpass  filter  such  as  a  a1-*1  order  Butterworth 
filter  could  be  used  to  eliminate  this  delay,  however  the  Graham  filter  was 
selected  for  the  results  presented  here  due  to  its  excellent  filter  properties 
and  availability.  The  1  to  2  second  delay  was  not  considered  significant  for 
the  real  time  monitoring  of  blade  modes. 

.The  filtered  output  y^  at  data  point  i  obtained  using  the  numerical 
convolution  is 


hj  yj+i 


(2.39) 


where  yj+£  is  the  raw  data  which  requires  M  data  points  ahead  and  M  data  points 
behind  point  i.  The  weight  function  hj  is  computed  from 


h, 

J 


u 

2jT 


(sin  o>tjT  +  sin  W^jT) 
it2  -  (u)t-coc)2(jT)2 


j  =  -M , . . . , +M 


j  *  0 


(2.40) 


and 


(J  +  w 

t  c 

2t t 


,  1  =  0 


The  cutoff  frequency  wc  and  termination  frequency  are  selected  to  achieve  the 
desired  passband  about  a  center  frequency  of  zero  (low  pass  filter).  The  weights 
hj  are  then  normalized  by  the  sum  of  all  the  weights.  To  achieve  bandpass  filtering 
about  a  center  frequency  u)0,  the  weights  are  shifted  by  multiplication  of  hj  by 
cos  UJ0jT.  Eq*  (2.39)  and  (2.40)  describe  the  low  pass  filter  (ref.  13)  and  the  bandpass 
filter  is  obtained  by  multiplication  of  the  weights  hj  by  cos  w0jT  and  is  described 
in  reference  14. 
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The  use  of  the  bandpass  filter  on  both  simulation  data  and  test  data  is 
presented  in  the  next  section  along  with  the  recursive  maximum  likelihood  iden¬ 
tification  results. 


t 


RESULTS  FROM  SIMULATION  ANALYSIS 

Before  analysis  of  test  data  is  performed  it  is  important  to  assess  the 
characteristics  and  limitations  of  the  identification  method  by  simulation.  Of 
particular  importance  is  the  convergence  characteristics,  (i.e.  the  quantifica¬ 
tion  of  parameter  accuracy  as  the  data  is  processed)  The  amplitude  of  the 
random  response,  the  initial  parameter  start-up,  and  the  digital  filter  are 
investigated  as  to  their  effect  on  identified  parameter  accuracy.  In  addition, 
the  standard  deviation  of  the  parameter  estimate  as  computed  in  the  algorithm 
must  be  evaluated  as  to  its  validity  as  a  measure  of  identified  parameter 
accuracy. 


Simulation  Model 

The  model  used  in  the  simulation  analysis  is  intended  to  represent  two  situa¬ 
tions.  First,  a  moderately  damped  mode  (t,  =  10%  critical  damping)  representative 
of  an  articulated  rotor  or  highly  damped  hingless  rotor  is  investigated.  Second, 
a  lightly  damped  mode  (£  =  2%  critical  damping)  representation  of  a  hingless  rotor 
is  investigated.  In  both  cases,  the  simulation  model  represents  a  single  uncoupled 
mode.  Such  a  mode  can  occur  as  an  isolated  blade  mode  or  a  lightly  damped  "system" 
mode  such  as  the  in-plane  regressing  mode  typical  of  ground  resonance.  Although 
the  identification  procedure  is  applicable  to  any  number  of  coupled  modes,  the 
results  presented  are  for  the  single  mode  case  which  is  perhaps  the  most  important 
for  rotor  stability  assessment.  In  addition,  the  modal  identification  procedure 
can  be  repeated  any  number  of  times  on  the  same  data,  in  each  case  a  different 
mode  isolated  with  the  bandpass  filter. 

The  model  used  in  the  simulation  is  a  state  space  second  order  system  and 
is  represented  by 


x^  =  x2  (3.1) 

2 

x0  =  -oj  x,  -  2^oo  x„  +  w(t)  (3.2) 

i  n  i  n  z 

where  x^  is  the  blade  modal  response  and  w(t)  is  zero-mean  white  gaussian  process 
noise. 

The  modal  measurement  is  obtained  from 

y  =  xl  +  r  (3.3) 

where,  r  is  zero-mean  white  gaussian  measurement  noise. 

The  stochastic  model  given  by  Eq.  (3.1)  through  Eq.  (3.3)  is  simulated  in 
discrete  form  by 

x1(i+l)  =  x^(i)  +  x2(i)  At  (3.4) 

x2(i+l)  =  — At  x^(i)  +  (l-2£wn  At)  x2(i)  +  wd(i) 
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and  y(i)  =  x^i)  +  r(i)  (3.5) 

where  is  zero-mean  white  discrete  process  noise  and  r  is  zero-mean  white  discrete 
measurement  noise. 

In  order  to  obtain  acceptable  integration  accuracy  using  the  discrete  model 
of  Eq.  (3.4)  and  Eq.  (3.5),  the  integration  step  size  At  should  be  at  least  10 
times  smaller  than  the  data  sample  time  between  the  i  and  i+1  sample.  For  all 
simulation  runs  At  was  selected  as  At  =  .002  seconds  and  the  measurement  y(i)  was 
made  every  At  =  .02  seconds. 

The  maximum  likelihood  algorithm  identifies  the  coefficients  of  the  equivalent 
ARMA  representation  as  discussed  in  the  previous  section.  The  ARMA  representation 
used  is 


y  +  a  y  +  a  y.  =  c  v  +  c  v.  + 
l  1  l-i  i  i — z  1  l-i  l  i-2 


v, 


(3.6) 


As  discussed  previously,  Eq.  (3.6)  is  a  single  input-single  output  discrete 
model  and  is  equivalent  to  an  innovations  or  Kalman  filter  representation  of  the 
system  described  by  Eq.  (3.4)  through  Eq.  (3.5).  The  filter  gains  are  embedded 
in  the  c  coefficients  of  Eq.  (3.6).  The  a  coefficients  are  the  coefficients  of 
the  characteristic  equation  from  which  the  damping,  £,  and  frequency,  tOjj,  can  be 
computed  as  shown  in  Eq.  (2.20)  through  Eq.  (2.25)  of  the  previous  section. 

Moderate  Damping  Case,  £  =  10%,  =  23.15  rad/sec 

The  simulation  model  shown  by  Eq.  (3.4)  and  Eq.  (3.5)  was  used  to  generate 
5000  data  samples.  The  measurement  noise  standard  deviation  used  is  ar  -  .1  and 
various  levels  of  process  noise  were  investigated.  The  data  used  in  the  identifi¬ 
cation  algorithm  is  unfiltered.  The  initial  parameter  start-up  is 

0(0)  =  {-1.5  .80  0.  0.)  (3.7) 


The  parameter  a£  =  .80  represents  a  critical  damping  of  greater  than  20%. 

The  initial  parameter  covariance  used  is 

P(0)  =  {.1  .1  1.0  1.0}  (3.8) 


Three  different  levels  of  process  noise  were  simulated  (o  =  1.0,  .2,  and  .04). 
This  corresponds  to  a  state  response  with  RMS  levels  of  1.25  deg,  .25  deg,  and  .05 
deg,  respectively.  Identified  parameter  convergence  is  shown  in  figure  2  through 
figure  4  for  the  three  levels  of  random  excitation.  As  shown  in  the  figures,  con¬ 
vergence  to  the  true  parameter  values  occurrs  in  each  case,  with  slower  convergence 
achieved  for  the  smaller  levels  of  process  noise  excitation.  The  ARMA  coefficients 
are  shown  on  the  left  hand  scale  and  the  damping  ratio  is  shown  on  the  right  hand 
scale  of  the  plot  for  coefficient  &2-  The  damping  ratio  is  computed  according  to 


G 


1 

2u)  At 

n 


Ml/a2) 


(3.9) 


The  identified  c  coefficients  account  for  the  presence  of  measurement  noise. 

If  the  ratio  of  process  to  measurement  noise  is  large  as  in  figure  2  (o^/cr  =  100), 

w  r 
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2  2 

Chen  the  c  coefficients  are  small.  As  the  ratio  o  /o  decreases,  the  c  coefficients 

w  r  2  2 

increase  in  magnitude  as  shown  in  figure  4.  In  the  limit  as  0  /a  approaches  °°  then 

w  r 

c  approaches  zero,  which  reduces  to  a  least  square  estimator  for  the  a  coefficients. 
Thus,  the  c  coefficients  in  the  maximum  likelihood  method  accounts  for  the  presence 
of  both  measurement  and  process  noise  by  automatically  "tuning"  the  Kalman  filter  in 
proportion  to  the  ratio  of  process  to  measurement  noise. 

Effect  of  Initial  Parameter  Estimate.  The  parameter  identif ication  results 
just  presented  for  the  10%  critical  damping  case  is  repeated  for  different  initial 
parameter  estimates.  The  second  set  for  the  initial  parameter  estimate  and  covariance 
is 


0(0)  =  (0  0  0  0}  (3.10) 

P(0)  =  {1.0  1.0  1.0  1.0}  (3.11) 


The  identified  parameter  convergence  is  shown  in  figure  5  and  figure  6  for 
process  noise  levels  of  Ow  =  .2  and  Ow  =  .04,  respectively.  Again,  identified 
parameter  convergence  is  excellent.  In  all  cases  studied  the  final  convergence 
was  found  to  be  independent  of  the  initial  parameter  start-up. 

Lightly  Damped  Case,  Z,  =  2%,  u)  =  23.15  rad/sec 

n 

The  simulation  was  repeated  using  a  critical  damping  of  Z,  =  2%.  This  case  is 
representative  of  a  lightly  damped  hingeless  rotor.  Two  levels  of  process  noise 
were  used  (Cw  =  .2  and  Ow  =  .04).  This  random  excitation  produces  a  RMS  random 
response  of  .75  deg  and  .15  deg,  respectively.  The  identified  parameter  convergence 
for  these  two  cases  is  shown  in  figure  7  and  figure  8.  Parameter  convergence  is 
excellent  and  occurs  more  rapidly  than  for  the  10%  damping  case  (compare  fig.  7 
with  fig.  3  and  fig.  8  with  fig.  4).  The  RMS  of  the  random  response  is  larger  for 
the  low  damping  case,  as  is  expected,  and  the  damping  is  more  rapidly  identified. 

Data  Length  Required  ^or  Accurate  Identification 

Based  upon  the  limited  number  of  simulation  runs  performed  an  estimate  of  the 
number  of  data  samples  required  for  accurate  identification  can  be  approximated.  As 
was  shown  in  the  simulation  study,  the  required  data  length  for  convergence  depends 
on  the  process  noise  level  and  damping  value.  Since  the  actual  process  noise  level 
is  not  known  in  a  wind  tunnel  test  program,  the  RMS  amplitude  of  the  random  response 
is  presented  rather  than  process  noise  level  (c^) .  Figure  9  presents  a  graphical 
summary  of  the  simulations  performed  and  shows  the  number  of  data  samples  required 
to  identify  damping  ratio  (%  critical)  to  within  10%  error.  This  is  plotted  vs  the 
RMS  amplitude  of  random  response  in  degrees.  Two  different  damping  ratios  (C  =  10% 
and  Z,  »  2%)  are  presented.  The  results  shown  in  figure  9  are  obtained  for  a  measure¬ 
ment  noise  level  of  <7r  »  .1  degrees,  u)n  =  23.15  rad/sec,  and  At  =  .02  seconds  sample 
time  interval. 

The  number  of  samples  required  is  very  dependent  upon  the  RMS  amplitude  and 
the  damping  ratio.  Typically,  the  number  of  samples  required  ranges  from  1000 
samples  for  low  damping  (2%)  and  large  RMS  amplitude  response  (1.25  deg^  to  greater 
than  5000  samples  for  high  damping  (10%)  and  low  RMS  response  (.05  degrees).  The 
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results  in  figure  9  should  be  used  only  as  a  guide  since  the  results  are  presented 
for  one  value  of  measurement  noise  (ar  =  .1),  using  unfiltered  data,  and  a  sampling 
rate  of  At  =  .02  seconds. 

The  results  of  figure  9  show  that  identified  parameter  convergence  requires 
between  20  seconds  to  100  seconds  of  data  (At  =  .02)  depending  upon  the  amplitude 
of  random  response  and  system  damping  level.  A  significant  finding  is  that  damping 
can  be  accurately  identified  from  random  response  levels  as  low  as  .05  degrees  in 
the  presence  of  .1  degrees  measurement  noise.  The  smaller  the  RMS  response,  the 
more  data  samples  required  for  parameter  convergence. 

Digital  Filter  Selection 

A  bandpass  filter  is  required  to  remove  unwanted  modes  and  1  per  rev  contamina¬ 
tions  from  the  data  before  identification  is  performed.  This  section  demonstrates 
the  effect  of  the  bandpass  filter  on  identification  accuracy  for  various  filter 
designs.  Figure  10  shows  the  amplitude  vs.  frequency  characteristics  of  four 
digital  filters  used  to  separate  the  desired  mode.  The  mode  of  interest,  uy,  occurs 
at  3.68  Hz  (23.15  rad/sec)  and  the  1  per  rev  (lp)  signal  is  shown  at  11.17  Hz. 

Filter  number  1,  2  and  4  of  figure  10  are  bandpass  filters,  whereas  filter  number  3 
is  a  low  pass  filter. 

The  simulation  data  was  generated  with  a  process  noise  level  of  aw  =  .2  for 
the  2%  critical  damping  case.  Four  separate  runs  were  made  and  identified  para¬ 
meter  convergence  was  assessed.  Figure  11  shows  identified  parameter  convergence 
for  filter  number  4.  The  identified  damping  value  is  somewhat  biased  from  the 
true  value  of  4  =  2%.  Filter  number  1,  2  and  3  showed  a  larger  bias,  indicating 
that  some  of  the  signal  (i.e.  the  random  response)  was  removed  by  the  filter.  The 
smaller  the  filter  bandwidth,  the  more  signal  removed  and  thus,  the  larger  the 
error.  Figure  12  shows  the  effect  of  filter  bandwidth  on  identified  damping  ratio 
for  the  four  filters  used. 

It  is  shown  that  a  filter  bandwidth  between  8  and  10  Hz  is  required  to  achieve 
acceptable  accuracy.  Normalizing  the  frequency  by  the  natural  frequency,  the  following 
approximate  frequency  range  should  be  selected  for  the  filter  assuming  a  white  noise 
random  excitation 


•27w  <w<  2.4w  (3.12) 

'-y  S 

The  relationship  in  Eq.  (3.12)  is  approximate  and  in  practice  the  filter  band¬ 
width  should  be  chosen  as  large  as  possible  within  the  0  to  Ip  frequency  band  for 
modes  below  1  p.  For  modes  above  lp  but  below  2p  the  filter  should  be  selected 
between  lp  and  2p  frequencies. 

For  a  mode  which  lies  very  close  to  the  lp  frequency,  it  may  be  required  to 
remove  the  lp  frequency  with  a  separate  estimator  before  using  the  bandpass 
filter. 

Figure  11  shows  that  the  identified  a  coefficients  are  somewhat  biased  by 
the  filter.  The  c  coefficients  are  significantly  biased.  This  is  because  the 
c  coefficients  are  influenced  by  the  ratio  of  process  noise  to  measurement  noise 
which  assumes  the  random  noise  is  white  noise.  Since  filtering  produces  colored 
innovations,  this  is  reflected  in  the  identified  c  coefficients. 
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Covariance  of  the  Parameter  Estimates 


Real  time  identification  of  modal  damping  and  frequency  is  not  complete  with¬ 
out  knowledge  of  the  uncertainty  in  the  paramr  er  estimate.  The  recursive  maximum 
likelihood  algorithm  provdies  a  measure  of  the  parameter  covariance  by  use  of  the 
Cramer-Rao  lower  bound.  A  discussion  of  theorems  relating  the  accuracy  of  maximum 
likelihood  estimate  to  the  Cramer-Rao  lower  bound  is  given  in  reference  12.  The 
Cramer-Rao  bound  is  a  lower  bound  on  estimation  accuracy  which  is  achieved  by  the 
maximum  likelihood  estimator  asymptotically.  An  estimate  of  the  parameter  covariance 
is  obtained  Dy 

cov(0N)  =  PN  (3.13) 


where  P^  is  obtained  from  Eq. 
from 


(2.35)  and  the  innovations  covariance  is  computed 


o  1  N 

xn-i£ 


i=l 


(3.14) 


Since,  in  practice,  the  noise  characteristics  of  the  innovations  is  often 
found  to  be  "colored"  rather  than  white,  the  parameter  covariance  as  given  by  Eq. 
(3.13)  must  be  corrected  for  non-white  innovations  (ref.  15).  This  can  be 
accomplished  using 

cov(0)  =  cov ( 0)  •  (y~^  At)  (3.15) 

- A 

where,  cov(0)  is  the  covariance  corrected  for  non-white  innovations  and  BW  is  the 
bandwidth  of  the  colored  noise  which  is  set  equal  to  the  digital  filter  bandwidth 
of  the  bandpass  filter  used  on  the  data.  The  correction  term  is  approximate  and 
is  found  to  work  well  in  practice  as  discussed  in  reference  15. 

Since  cov  (0)  represents  the  covariance  of  the  coefficients  of  the  charac¬ 
teristic  equation,  a  transformation  is  required  to  obtain  the  covariance  of  the 
damping.  Since  the  damping  ratio  is  related  to  the  coefficient  a2  by 

g  =  e2^ht  (3.16) 

the  covariance  of  £  can  be  obtained  by  lineariiation  of  Eq.  (3.16)  to  first  order 
as 


2,  *  1  -  2Cw  At 

2  n 

/A 

The  covariance  of  £  can  then  be  computed  from  Eq.  (3.17)  or 

A  _  A  2 
cov(0  =  cov(a„)/(2  At  u  ) 

2  n 

A 

The  standard  deviation  of  Z,  is  thus 

» 'csva'j 

op  »  - 

4  2  At  u 


(3.17) 


(3.18) 


(3.19) 
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Eq.  (3.15)  is  used  to  compute  the  standard  deviation  of  the  parameters  in 
the  ARMA  model  (e.g.  the  a-  coefficients)  and  Eq.  (3.19)  is  used  to  estimate  the 
standard  deviation  of  the  estimated  damping.  Figure  13  shows  the  identified  para¬ 
meter  convergence  and  standard  deviation  of  the  estimate  for  a^  and  &2  f°r  a 
simulation  damping  ratio  of  2%  and  aw  =  .2.  The  estimated  standard  deviation  pro¬ 
vides  an  excellent  measure  of  uncertainties  for  the  identified  parameters. 

Figure  14  shows  identified  parameter  convergence  and  +  2a  bands  for  a  simula¬ 
tion  damping  of  2%  with  process  noise  0^  =  .04.  The  la  bands  were  slightly 
optimistic  while  the  2a  bands  always  contains  the  true  damping  value.  Conver¬ 
gence  is  slower  for  figure  14  than  figure  13  due  to  the  smaller  RMS  of  the  random 
response  (aw  =  . 04  vs  aw  =  . 2) .  In  all  simulation  studies  performed  the  +  la 
band  or  +  2a  band  provided  an  accurate  measure  of  the  parameter  uncertainty. 


RESULTS  FROM  TEST  DATA 


Test  Data  Description 

A  limited  amount  of  test  data  obtained  from  a  4-bladed  hingeless  rotor  test 
in  the  NASA  Langley  Research  Center  (LRC)  Helicopter  Hover  Facility  in  1981  was 
furnished  by  Mr.  William  T.  Yeager,  Jr.  and  Mr.  Wayne  R.  Mantay.  This  data  is 
used  to  identify  in-plane  mode  damping  and  frequency.  The  rotor  was  operated  in 
hover  on  the  Aeroelastic  Rotor  Experimental  System  (ARES)  over  a  rotor  speed  range 
of  0  to  618  RPM.  The  Helicopter  Hover  Facility  is  used  for  hover  testing  on  the 
ARES  before  testing  in  the  LRC  Transonic  Dynamics  Wind  Tunnel.  Two  RPM  condi¬ 
tions,  250  RPM  and  618  RPM,  were  investigated. 

Of  particular  concern  during  rotor  speed  increase  is  the  potential  for  ground 
resonance  which  occurs  near  the  rotor  in-plane  regressing  mode  frequency  (Q-ojr ) . 
Since  the  in-plane  lead-lag  angle  measurement  occurs  in  the  rotating  system  the 
bandpass  filter  was  chosen  to  be  centered  around  the  rotating  in-plane  frequency 
U3J-.  The  250RPM  condition  was  observed  to  have  a  large  resonance  close  to  the 
rotor  RPM,  excited  purely  by  random  turbulence.  The  in-plane  magnitude  response 
was  approximately  +  1  degree.  Because  of  this  large  resonance  the  damping  was 
expected  to  be  quite  small  and  a  forced  control  input  excitation  would  be  dangerous 
to  apply.  The  recursive  maximum  likelihood  approach  to  identify  modal  damping  from 
random  response  is  thus  ideally  suited  for  this  RPM  condition. 

The  618  RPM  condition  was  considerably  more  damped  and  the  modal  response 
due  to  inherent  random  excitation  was  less  than  .05  degrees.  It  is  expected 
that  a  longer  data  record  length  would  be  required  for  the  618  RPM  condition 
for  accurate  identification  because  of  the  small  amplitude  response. 

The  test  data  used  for  the  250  RPM  and  618  RPM  conditions  are  shown  in  figure 
15.  Figure  15  shows  only  a  portion  of  the  50  seconds  of  test  data  record  length 
used  in  the  identification.  The  lead-lag  measurement  for  the  250  RPM  case  of 
figure  15  shows  the  +  1  degree  resonance  response  occurring  at  in-plane  frequency 
wr.  The  sinusoidal  response  shown  in  figure  15  (bottom  figure)  for  the  measure¬ 
ment  of  lead-lag  at  the  618  RPM  condition  is  due  to  the  Ip  frequency  of  the  rotor 
and  is  not  the  in-plane  mode.  The  large  lp  response  occurs  at  this  condition  due 
to  a  non-zero  trim  collective  and  cyclic  control  (0O  =  8°,  Bq  =  4°).  The  in-plane 
mode  amplitude  is  less  than  .05  degrees  which  will  be  discussed  in  more  detail  later. 

Identification  Results  at  250  RPM 

The  data  at  the  250  RPM  condition  shown  in  figure  15  (top  figure)  was  processed 
with  bandpass  filter  number  4,  the  characteristics  of  which  are  shown  in  figure  10. 
During  the  process  of  transforming  the  analog  data  tape  to  digital  format  a  portion 
of  the  data  was  distroyed.  The  region  of  bad  data  occurred  between  the  1200th  and 
1600th  samples.  The  complete  data  length  processed  consisted  of  2500  samples  at  .02 
seconds  per  sample  rate. 

The  filtered  data  was  used  in  the  recursive  maximum  likelihood  algorithm  and 
the  ARMA  model  coefficients  identified.  The  results  are  shown  in  figure  16  for 
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the  coefficients  of  the  characteristic  equation,  a^  and  £2  and  the  damping  and 
frequency  estimates  are  shown  in  figure  17,  The  region  of  bad  data  shown  in 
these  figures  was  also  recursively  processed  by  the  algorithm.  Identified 
parameter  convergence  is  more  accurate  before  the  bad  data  region  is  encountered 
(i.e.  up  to  1200  samples). 

Convergence  is  shown  to  be  excellent  in  figure  16  and  17  as  indicated  by  the 
+  2a  accuracy  bands  shown  in  the  figure  before  the  region  of  bad  data  is  encountered. 
The  identified  damping  ratio  at  1200  samples  is  1.4%  with  1  standard  deviation  equal 
+  *7%.  The  rapid  convergence  rate  is  a  result  of  the  large  amplitude  random  re¬ 
sponse  (+  1  degree)  as  observed  in  the  test  data.  The  convergence  rate  shown  in 
figure  16  and  17  is  consistent  with  that  expected  from  the  simulation  studies. 

For  example,  referring  to  figure  9,  it  is  shown  that  perhaps  1500  data  points  are 
required  to  achieve  the  desired  accuracy.  Ideally,  it  would  be  desirable  to  pro¬ 
cess  the  full  2500  samples  of  data,  however,  because  of  the  bad  data  region,  the 
parameter  estimate  was  shifted  as  shown  in  figure  16  and  17.  Notice  however, 
even  with  bad  data,  the  algorithm  is  performing  properly,  in  that  the  parameter 
estimates  are  attempting  to  recover  as  good  data  becomes  available. 

The  identified  frequency  is  shown  in  figure  17  and  the  estimated  value  at 
1200  samples  is  26.3  rad/sec.  This  value  is  nearly  the  same  as  determined  by 
that  obtained  from  shake  test  data  and  analytic  predictions. 

The  identified  damping  and  frequency  parameters  shown  in  figure  17  demonstrates 
that  the  identification  procedure  is  performing  successfully,  and  is  in  excellent 
agreement  to  that  obtained  in  the  simulation  studies. 

Identification  Results  at  618  RPM 

The  hover  test  data  obtained  at  the  618  RPM  condition  is  shown  again  in 
figure  18  with  an  expanded  scale.  Figure  18  shows  the  unfiltered  data  (top 
figure)  and  the  bandpass  filtered  data  (lower  figure).  The  unfiltered  data 
shows  the  dominant  response  is  the  lp  rotor  frequency  (+  .3  degrees).  The 
filtered  data  removes  the  lp  frequency  and  retains  the  desired  inplane  mode 
which  occurs  at  u)j-  =  38.83  rad/sec.  This  mode  is  clearly  shown  in  the  lower 
figure  of  figure  18.  The  frequency  of  this  mode  appears  to  change  frequency 
throughout  the  data  length  and  is  due  to  the  random  excitation  driving  the  rotor 
system.  The  RMS  amplitude  of  response  is  approximately  .03  degrees  as  shown  in 
the  lower  figure.  Referring  back  to  figure  10,  the  simulation  study  has  shown 
that  5000  data  points  at  a  sample  rate  of  .02  seconds  (100  seconds  of  data)  are 
required  for  acceptable  accuracy  for  a  RMS  random  response  of  .05  degrees. 

Since  only  45  seconds  of  data  are  available  at  the  618  RPM  condition  and  the  RMS 
amplitude  is  .03  it  is  expected  that  convergence  would  require  approximately  twice 
as  much  data.  This  was  in  fact  the  case  and  the  convergence  plots  are  shown  in 
figure  19.  The  identified  damping  ratio'  r,,  is  shown  in  the  lower  figure  and  the 
+  2a  accuracy  band  shows  a  longer  data  record  is  required  for  convergence.  The 
identified  damping  is  estimated  after  45  seconds  of  data  (4500  samples  at  .01 
seconds  sample  rate)  to  be  $  =  9%  with  1  standard  deviation  of  +  3.5%.  Indepen¬ 
dent  estimates  of  the  damping  at  this  condition  using  the  log  decrement  approach 
(visual  inspection  of  an  initial  condition  decay  response  was  used)  during  the 
hover  tests  reveal  a  damping  estimate  of  approximately  4%.  The  results  of  figure 
19  thus  contain  this  value  within  the  accuracy  band  (+  2o)  as  shown  in  this  figure. 
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The  identified  frequency  convergence  is  also  shown  in  figure  19  (top  figure) 
and  the  estimated  frequency  is  ui  =  38.83  radians/sec.  This  value  is  close  to 
that  obtained  from  independent  shake  test  data  and  analytic  predictions.  The 
frequency  parameter  is  identified  more  rapidly  than  the  damping  and  requires 
less  data  length.  This  was  found  to  be  the  case  in  the  simulation  study  per¬ 
formed  prior  to  analyzing  the  test  data. 

The  identification  results  at  the  618  RPM  condition  have  correctly  identified 
the  frequency  and  approximately  twice  the  data  length  available  is  required  to 
identify  the  damping  with  acceptable  accuracy.  This  conclusion  was  predicted  by 
the  simulation  studies  which  show  that  100  seconds  of  data  is  needed  (twice  the 
data  length  available)  to  obtain  the  required  accuracy  of  the  damping  estimate. 

Table  1  provides  a  summary  of  the  identified  damping  and  frequency  at  both 
the  250  and  618  RPM  conditions.  Comparisons  are  made  with  other  available  pre¬ 
diction  methods.  Shake  test  data  and  analytic  simulation  is  referred  to  in  the 
table  as  Analytic  Method.  The  Log  Decrement  method  is  based  upon  visual  inspec¬ 
tion  of  the  amplitude  decay  resulting  from  removal  of  a  forced  excitation.  The 
damping  values  obtained  from  the  Analytic  method  are  estimates  based  upon  the 
estimated  rotor  blade  damping.  The  actual  rotor  blade  damping  value  was  not 
precisely  known  (1%  critical  damping  was  used  in  the  simulation).  The  fixed 
system  damping  and  frequency  are  for  the  in-plane  regressing  mode  and  the 
rotating  system  values  are  shown  for  the  maximum  likelihood  method  and  log  decre¬ 
ment  method.  The  fixed  system  values  are  shown  for  the  analytic  method  and 
rotating  system  values  are  obtained  via  transformation  as  noted  in  the  table. 

From  the  data  available  at  the  250  RPM  and  618  RPM  condition  the  recursive 
maximum  likelihood  approach  performed  as  expected  from  the  simulation  studies. 

Identified  damping  convergence  to  acceptable  accuracy  depends  upon  the 
magnitude  of  the  RMS  random  response  and  the  level  of  modal  damping.  The  larger 
the  random  response  the  faster  convergence.  Typically  from  30  to  100  seconds  of 
data  are  required  for  acceptable  accuracy.  For  +  1  degree  random  response  and 
low  damping  (1.4%)  30  to  50  seconds  of  data  are  required  (250  RPM  condition). 

For  +  .03  degree  random  response  and  damping  approximately  4%,  approximately 
100  seconds  of  data  are  required  (618  RPM  condition). 
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CONCLUSIONS 


A  recursive  maximum  likelihood  identification  procedure  has  been  presented 
for  the  estimation  of  helicopter  rotor  blade  modal  damping  and  frequency  from 
random  response  data.  The  algorithm  was  shown  to  have  excellent  convergence 
properties,  provide  accurate  estimates  of  parameter  covariance,  and  require 
virtually  no  user  interaction.  As  a  result,  the  RML  method  can  be  used  for 
continuous  estimation  of  damping  and  frequency  during  wind  tunnel  (or  flight) 
testing. 

Simulation  studies  have  demonstrated  that  convergence  to  true  parameter 
estimates  occur  even  for  very  small  amplitude  RMS  random  response.  RMS  random 
response  amplitude  of  .05  degrees  to  1.25  degrees,  wnich  are  typical  of  that 
expected  in  wind  tunnel  data,  resulted  in  successful  identification.  Data  record 
lengths  required  for  accurate  identification  range  from  20  seconds  to  100  seconds 
of  data  depending  upon  RMS  amplitude  of  random  response  and  damping  level.  These 
data  lengths  are  compatible  with  typical  wind  tunnel  test  data  taken  at  steady 
rotor  trim  conditions.  Since  a  forced  input  is  not  required,  damping  and  frequency 
can  be  identified  without  interfering  with  other  wind  tunnel  test  objectives.  The 
estimated  parameter  covariance  was  shown  to  provide  accurate  confidence  hands  on 
the  damping,  thus  providing  a  indicator  as  to  length  of  data  required  at  a  given 
test  condition.  Since  identification  parameter  convergence  was  shown  to  be 
independent  of  initial  parameter  start-up  values,  the  RML  procedure  can  be  used 
with  virtually  no  user  interaction.  This  proved  to  be  the  case  on  all  data  pro¬ 
cessed  in  this  investigation. 

The  hover  facility  test  data  identified  parameter  convergence  results  obtained 
were  consistent  with  the  simulation  studies.  At  250  RPM,  identified  in-plane  mode 
parameters  are  £  =  1.4%  critical  damping  and  =  26.3  rad/sec.  The  +  lo  estimate 
on  r,  is  .7%  using  24  seconds  of  data.  The  identified  frequency  is  in  agreement 
with  independent  shake  test  data  results.  Accurate  estimates  of  damping  from  other 
independent  tests  were  not  available.  At  618  RPM,  the  identified  frequency  was 
u >r  =  38.83  rad/sec  which  is  in  agreement  with  independent  shake  test  data.  The 
identified  damping  is  estimated  using  45  seconds  of  data  to  be  r.  =  9%  with 
+  la  =  3.5%.  The  large  standard  deviation  shows  that  a  longer  data  record  is 
required  which  is  consistent  with  the  simulation  results  which  showed  that  twice 
the  data  length  is  needed  to  obtain  an  accurate  estimate.  For  both  the  250  RPM 
and  618  RPM  the  RML  identified  parameter  convergence  performed  as  predicted  by 
the  simulation  studies. 


RECOMMENDATIONS 


The  emphasis  of  this  investigation  was  on  a  single  mode  response  separated 
from  other  modes  with  a  bandpass  filter.  Although  the  techniques  developed  are 
applicable  to  any  number  of  over-lapping  modes,  simulation  studies  on  test  data 
analysis  have  only  been  performed  on  a  single  mode  case.  Since  rotorcraft  rotor 
instability  can  result  from  coupled  modes  it  is  recommended  that  simulation 
studies  analogous  to  that  performed  in  this  report  be  repeated  for  coupled  mode 
response  of  two  and  possibly  three  modes.  In  addition,  criteria  should  be 
investigated  to  automatically  determine  the  required  number  of  modes  in  the  data. 

The  RML  technique  contains  provision  to  track  time  varying  parameters. 

Simulations  should  be  performed  to  quantify  this  capability.  This  would  be 
useful  for  monitoring  damping  changes  while  going  from  one  steady  rotor  trim 
condition  to  another  trim  condition.  The  RML  technique  contains  provision  for  identi¬ 
fication  of  data  generated  by  known  forced  control  inputs  and  free  response  data. 

Both  these  situations  have  not  been  investigated  via  simulation.  Simulation  studies 
should  be  performed  to  quantify  this  capability. 

Finally,  further  test  data  applications  should  be  performed.  Both  on-line  and 
off-line  applications  would  provide  a  data  base  for  further  assessment  of  the  RML 
technique. 
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Table  1.  -  Test  Data  Identified  Damping  and  Frequency  and 
Comparisons  with  Analytic  and  Log-Decrement 
Predictions 


Estimation 

Rotating  System 

Fixed 

System 

Method 

GJ 

r 

SR 

Rad/Sec 

?R 

%  Critical 

^F 

Rad /Sec 

?F 

%  Critical 

n  =  250  RPM 

Analytic  Method ^ 

26.2(3) 

0.55%(A) 

.8 

18% 

Log  Decrement 

Not 

Used^- 

— 

— 

Max.  Likelihood 

26.3 

1.4%(la=. 7%) 

— 

— 

Q  =  618  RPM 

Analytic  Method ^ 
(2) 

Log  Decrement 

37.63(3) 

0.33 %(4) 

4% 

25.9 

.5% 

Max.  Likelihood 

38.83 

9% (la-3. 5%) 

— 

— 

(1)  Based  upon  analytic  simulation  and  shake  test  data. 

(2)  Computed  based  upon  visual  inspection  of  amplitude  decay  response 
resulting  from  removal  of  forced  excitation. 


(3) 

(4) 

(5) 


Computed  from  shake  test  data. 

0) 


Computed  from  C 

K 


Forced  excitation  riot  used  at  250  RPM  because  of  large  observed 
resonance  at  this  condition. 


(6)  Values  shown  are  obtained  at  slightly  different  RPM. 
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Test  Data 


Figure  1. 


C,  W  ,  and 
Standard  Deviations 


Filter  Desired 
Modal  Response 


Identify  Coefficients 
of  Characteristic 
Equation 


Calculate  Z-Transform 
Domain  Roots 


Calculate  Continuous 
Time  Damping 
and  Frequency 


-  Overall  Procedure  for  Real  Time  Identification 
of  Modal  Damping  and  Frequency. 
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£  £  =  10%  Critical  Damping 
0  ^  =  2%  Critical  Damping 


RMS  Amplitude  of 


Random  Response,  deg 


Figure  9.  * 


•  Pd  To  Identify  Damping  Ratio  to 
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Figure  11.  -  Identified  Parameter  Convergence  Using  Digital  Filter 
Number  4  (RMS  of  Random  Response  =  .75  deg,  a^.l  deg, 
O  “.2  deg/sec^,  At*. 02  sec,  u>  a23.l5  rad/sec,  C=2%) . 
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Figure  12.  -  Effect  of  Filter  Bandwidth  on  Identified  Damping  Ratio 

(RMS  of  Random  Response  *  .75  deg,  ar=.l  deg,  aw=.2  deg/sec  , 
At-. 02  sec,  u)  =23.15  rad/sec,  £-27,) . 
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Figure  16.  -  Identified  Parameter  Convergence  From  Hover  Facility  Test 

Data  at  fi  =  250  RPM  (Parameter  Accuracy  is  Best  before  Region 
of  Bad  Data  is  encountered  at  1200  samples,  At*. 02  sec). 
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Hover  Facility  Test  Data  at  ft  =  618  RPM  Showing  Unfiltered  and  Filtered  Test  Data 
(Only  a  portion  of  Data  is  shown,  Amplitude  Scale  is  Approximate  ,  0  =  8°,  B  =4°) 
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